Causal dissipative hydrodynamics for QGP fluid in 2+1 dimensions 
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In 2nd order causal dissipative theory, space-time evolution of QGP fluid is studied in 2+1 di- 
mensions. Relaxation equations for shear stress tensors are solved simultaneously with the energy- 
momentum conservation equations. Comparison of evolution of ideal and viscous QGP fluid, ini- 
tialized under the same conditions, e.g. same equilibration time, energy density and velocity profile, 
indicate that in a viscous dynamics, energy density or temperature of the fluid evolve slowly, than 
in an ideal fluid. Cooling gets slower as viscosity increases. Transverse expansion also increases in 
a viscous dynamics. For the first time we have also studied elliptic flow of 'quarks' in causal viscous 
dynamics. It is shown that elliptic flow of quarks saturates due to non-equilibrium correction to 
equilibrium distribution function, and can not be mimicked by an ideal hydrodynamics. 

PACS numbers: 47.75.+i, 25.75.-q, 25.75.Ld 



I. INTRODUCTION 

One of the most important discoveries in Relativistic 
Heavy ion collider (RHIC) at Brokhaven National Lab- 
oratory is the large elliptic flow in non-central Au-|-Au 
collisions [ESSH • Elliptic flow measures the momen- 
tum anisotropy of produced particles and is quantified by 
the 2nd harmonic of the azimuthal distribution, 



V2{pt) =< cos{2(j)) >= 



dN 



(1.1) 



dycPpT 



Elliptic flow is naturally explained in hydrodynamics. 
Hydrodynamic pressure is built up from rescattering of 
secondaries, and pressure gradients drive the subsequent 
collective motion. In non-central Au-|-Au collisions, ini- 
tially, the reaction zone is asymmetric (almond shaped). 
The pressure gradient is large in one direction and small 
in the other. The asymmetric pressure gradients gener- 
ates the elliptic flow. Naturally, in a central collision, 
reaction zone is symmetric and elliptic flow vanishes. 
Observed elliptic flow then give the strongest indication 
that in non-central Au-t-Au collisions, a collective QCD 
matter is produced. Whether the formed matter can be 
identified as the much sought after Quark-Gluon Plasma 
(QGP) as predicted in Lattice QCD simulations [5] is 
presently debatable. 

Ideal hydrodynamics has been partly successful in ex- 
plaining the observed elliptic flow, quantitatively ^ . El- 
liptic flow of identified particles, up to pt ~1.5 GeV are 
well reproduced in ideal hydrodynamics. Ideal hydrody- 
namics also explains the transverse momentum spectra 
of identified particles (up to pT ~ 1-5 GeV). Success of 
ideal hydrodynamics in explaining bulk of the data Q, 
together with the string theory motivated lower limit of 
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shear viscosity 77/s > l/47r to a paradigm 

that in Au-f Au collisions, a nearly perfect fluid is created. 

However, the paradigm of "perfect fluid" produced in 
Au-|-Au collisions at RHIC need to be clarifled. As in- 
dicated above, the ideal hydrodynamics is only partially 
successful and in a limited pT range [px <1.5 GeV) 0]. 
The transverse momentum spectra of identifled particles 
also starts to deviate form ideal fluid dynamics prediction 
beyond Pt ~ 1.5 GeV. Experimentally determined HBT 
radii are not reproduced in the ideal fluid dynamic mod- 
els, the famous "HBT puzzle" [l^]. It also do not repro- 
duce the experimental trend that elliptic flow saturates 
at large transverse momentum. These shortcomings of 
ideal fluid dynamics indicate greater importance of dis- 
sipative effects in the pr ranges greater than 1.5 GeV 
or in more peripheral collisions. Indeed, ideal fluid is a 
concept, never realized in nature. As suggested in string 
theory motivated models 0, Q , QGP viscosity could be 
small, rj/s > l/47r, nevertheless it is non-zero. It is im- 
portant to study the effect of viscosity, even if small, on 
space-time evolution of QGP fluid and quantify its effect. 
This requires a numerical implementation of relativistic 
dissipative fluid dynamics. Furthermore, if QGP fluid 
is formed in heavy ion collisions, it has to be charac- 
terized by measuring its transport coefficients, e.g. heat 
conductivity, bulk and shear viscosity. Theoretically, it 
is possible to obtain those transport coefficients in a ki- 
netic theory model. However, in the present status of 
theory, the goal can not be achieved immediately, even 
more so for a strongly interacting QGP (sQGP). Alter- 
natively, one can use the experimental data to obtain 
a " phenomenological" limit of transport coefficients of 
sQGP. It will also require a numerical implementation of 
relativistic dissipative fluid dynamics. There is another 
incentive to study dissipative hydrodynamics. Ideal hy- 
drodynamics depends on the assumption of local equilib- 
rium. Before local equilibrium is attained, the system has 
to pass through a non-equilibrium stage, where (if non- 
equilibrium effects are small) dissipative hydrodynamics 
may be applicable. Indeed, we can explore early times of 
fluid evolution better in a dissipative hydrodynamics. 
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Theory of dissipative relativistic fluid has been formu- 
lated quite early. The original dissipative relativistic fluid 
equations were given by Eckart and Landau and Lif- 
shitz 1^. They are called 1st order theories. Formally, 
relativistic dissipative hydrodynamics are obtained from 
an expansion of entropy 4-currcnt, in terms of dissipative 
fluxes. In 1st order theories, entropy 4-current contains 
terms linear in dissipative quantities. 1st order theory 
of dissipative hydrodynamics suffer from the problem of 
causality violation. Signal can travel faster than light. 
Causality violation is unwarranted in any theory, even 
more in a relativistic theory. The problem of causal- 
ity violation is removed in the Israel-Stewart's 2nd order 
theory of dissipative fluid [l3l ■ In 2nd order theory, ex- 
pansion of entropy 4-current contains terms 2nd order in 
dissipative fluxes. However, these leads to complications 
that dissipative fluxes are no longer function of the state 
variables only. They become dynamic. The space of ther- 
modynamic variables has to be extended to include the 
dissipative fluxes (e.g. heat conductivity, bulk and shear 
viscosity) . 

Even though 2nd order theory was formulated some 
30 years back, significant progress towards its numer- 
ical implementation has only been made very recently 
[11, M, M, M, M, [H, M, HH, m, m . At the Cyclotron 
Centre, Kolkata, we have developed a code " AZHYDRO- 
KOLKATA" to simulate the hydrodynamic evolution of 
dissipative QGP fluid. Presently only dissipative ef- 
fect included is the shear viscosity. Some results of 
AZHYDRO-KOLKATA, for first order dissipative hydro- 
dynamics have been published earlier [H, |20, HH . In the 
present paper, for the first time, we will present some 
results for 2nd order dissipative hydrodynamics in 2-1-1 
dimensions. In the present paper, we will consider effect 
of dissipation in the QGP phase only. Effect of phase 
transition will be studied in a later publication. 

The paper is organized as follows: In section [TT] we 
briefly review relativistic dissipative fluid dynamics. In 
section IIIII we derive the relevant equations in 2+1 di- 
mension (assuming boost-invariance) . Required inputs 
e.g. the equation of state, viscosity coefficient and initial 
conditions are discussed in section IIVI Simulation re- 
sults from AZHYDRO-KOLKATA are shown in section 
IVl In section rVIII we compare the transverse momentum 
spectra and elliptic flow of quarks in ideal and viscous 
dynamics. The concluding section IIXI summarizes our 
results. 



II. DISSIPATIVE FLUID DYNAMICS 

In this section, I briefly discuss the phenomenological 
theory of dissipative hydrodynamics. More detailed ex- 
position can be found in [l3| . 

A simple fluid, in an arbitrary state, is fully speci- 
fied by primary variables: particle current (N^), energy- 
momentum tensor (T^") and entropy current {S'^) and 
a number of additional (unknown) variables. Primary 



variables satisfies the conservation laws; 

d^N^ = 0, (2.1) 
d^Tf'" - 0, (2.2) 

and the 2nd law of thermodynamics, 

9^5"^ > 0. (2.3) 

In relativistic fluid dynamics, one defines a time-like 
hydrodynamic 4- velocity, (normalized as — 1). 
One also define a projector, A^'' — g^"^ — u'^u'^, orthog- 
onal to the 4- velocity (A'^''it,y = 0). In equilibrium, an 
unique 4-velocity (u^) exists such that the particle den- 
sity (n), energy density (e) and the entropy density (s) 
can be obtained from, 

K,^ nu, (2.4) 
= eu^u" - p^^" (2.5) 
^e^, = (2.6) 

An equilibrium state is assumed to be fully specified 
by 5-parameters, (n, e, u^) or equivalently by the thermal 
potential, a = fi/T (p being the chemical potential) and 
inverse 4-temperature, = u'^/T. Given a equation of 
state, s — s(e,n), pressure p can be obtained from the 
generalized thermodynamic relation, 

Sl!^ = pp^ ~ aNlf^ + f3^T,\^ (2.7) 

Using the Gibbs-Duhem relation, d{pP^) = Nj^^da — 
T^^dPx, following relations can be established on the 
equilibrium hyper-surface ^eqice, P^), 

dS^^ = ~adNlf^ + PxdT^^^ (2.8) 

In a non-equilibrium system, no 4-velocity can be 
found such that Eqs l2.4l2.5l2.6l remain valid. Tensor de- 
composition leads to additional terms, 

m = + SN^" = nu^" + V^" (2.9) 

= [eu^'u'' - pA""] + HA^" -f- nf"" 

+ {W''u'' + Wuf) (2.10) 
= Si^q -f 5Sf = su^' + (2.11) 

The new terms describe a net flow of charge = 
Af'^N^, heat flow, W = {e+p)/nV'' + (where 
is the heat flow vector), and entropy flow H = 

— ^Af^^T^" — p is the bulk viscous pressure and n^" = 
[i(A'^'"A''^ + A'^'^A'^^ - i A^''A'"^]T^^ is the shear stress 
tensor. Hydrodynamic 4-velocity can be chosen to elimi- 
nate either (the Eckart frame, is parallel to particle 
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flow) or the heat flow q^^ (the Landau frame, is par- 
allel to energy flow). In relativistic heavy ion collisions, 
central rapidity region is nearly baryon free and Lan- 
dau's frame is more appropriate than the Eckart's frame. 
Dissipative flows are transverse to and additionally, 
shear stress tensor is traceless. Thus a non-equilibrium 
state require 1-1-3-1-5=9 additional quantities, the dissi- 
pative flows n, q^^ (or V^) and n^'^ . In kinetic theory, 
and T'^" are the 1st and 2nd moment of the distri- 
bution function. Unless the function is known a-priori, 
two moments do not furnish enough information to enu- 
merate the microscopic states required to determine S''^, 
and in an arbitrary non-equilibrium state, no relation 
exists between, iV, T'"' and S"^. Only in a state, close 
to a equilibrium one, such a relation can be established. 
Assuming that the equilibrium relation Ea l2.8l remains 
valid in a "near equilibrium state" also, the entropy cur- 
rent can be generalized as, 



so in relativistic theory. Causality violation of dissipative 
hydrodynamics is corrected in 2nd order theories |il3i | . In 
2nd order theories, entropy current contain terms up to 
2nd order in the deviations, 7^ 0. The most general 
containing terms up to 2nd order in deviations can 
be written as, 

Q = -(Pon -Piq qu+P2T^v\T^ )^ 7^ h 

, (2.18) 

As before, one can cast the entropy production rate 
(Tdfj^S^^) in the form of Eq l2.14l Neglecting the terms 
involving dissipative flows with gradients of equilib- 
rium thermodynamic quantities (both are assumed to be 
small) and demanding that a linear relation exists be- 
tween the dissipative flows and thermodynamic forces, 
following relaxation equations for the dissipative flows 
can be obtained, 



5^ = S^g + rfS*^ = - aN" + f3xT^'' + Q" (2.12) 



where Q'^ is an undetermined quantity in 2nd order in 
deviations, SN^" = - N^^ and ST'''' = Tt'" - T^^ . 
Detail form of is constrained by the 2nd law d^S^' > 
0. With the help of conservation laws and Gibbs-Duhem 
relation, entropy production rate can be written as. 



d^S^ - -SNf'd^a + ST^^d^p, + d^Q^ (2.13) 

Choice of leads to 1st order or 2nd order theories 
of dissipative hydrodynamics. In 1st order theories the 
simplest choice is made, — 0, entropy current con- 
tains terms up to 1st order in deviations, 6N^ and ST^" . 
Entropy production rate can be written as. 



(2.14) 



where, X = -V.u; = ^ - u^d^u^' and X'"" = 

The 2nd law, 9^5^ > can be satisfied by postulat- 
ing a linear relation between the dissipative flows and 
thermodynamic forces, 



n = 

nv _ 



-A 



2?7 [V 



P T 



(2.19) 
(2.20) 
(2.21) 



where D = u'^d^is the convective time derivative. Unlike 
in the 1st order theories, in 2nd order theories, dynamical 
equations control the dissipative flows. Even if thermo- 
dynamic forces vanish, dissipative flows do not vanish 
instantly. 

Before we proceed further, it may be mentioned that 
the parameters, a and P\ are not connected to the actual 
state {N>',T'"'). The pressure p in Ea l^T^ is also not 
the "actual" thermodynamics pressure, i.e. not the work 
done in an isentropic expansion. Chemical potential a 
and 4-inverse temperature (3\ has meaning only for the 
equilibrium state. Their meaning need not be extended 
to non-equilibrium states also. However, it is possible to 
fit a fictitious "local equilibrium" state, point by point, 
such that pressure p in Ea l2.12l can be identified with 
the thermodynamic pressure, at least up to 1st order. 
The conditions of fit fixes the underlying non-equilibrium 
phase-space distribution. 



n= -c^, (2.15) 

q^= -A^V''(/i/r), (2.16) 

e +p 

tt"" = 2f]W<''u''> (2.17) 

where ^, A and 77 are the positive transport coefhcients, 
bulk viscosity, heat conductivity and shear viscosity re- 
spectively. 

In 1st order theories, causality is violated. If, in a given 
fluid cell, at a certain time, thermodynamic forces vanish, 
corresponding dissipative fluxes also vanish instantly. Vi- 
olation of causality is unwanted in any theory, even more 



III. (2+l)-DIMENSIONAL VISCOUS 
HYDRODYNAMICS WITH LONGITUDINAL 
BOOST INVARIANCE 

Complete dissipative hydrodynamics is a numerically 
challenging problem. It requires simultaneous solution of 
14 partial differential equations (5 conservation equations 
and 9 relaxation equations for dissipative flows). We re- 
duce the problem to solution of 6 partial differential equa- 
tions (3 conservation equations and 3 relaxation equa- 
tions). In the following, we will study boost-invariant 
evolution of baryon free QGP fluid, including the dissi- 
pative effect due to shear viscosity only. Shear viscosity 
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is the most important dissipative effect. For example, in 
a baryon free QGP, heat conduction is zero and we can 
disregard Eq l2.20l Bulk viscosity is also zero for the QGP 
fluid (point particles) and Ea l2.19l can also be neglected. 
Shear pressure tensor has 5 independent components but 
the assumption of boost invariance reduces the number 
of independent components to three. For a baryon free 
fluid, we can also disregard the conservation equation 
Ea l2.1l With the assumption of boost-invariance, energy- 
momentum conservation equation d^T^"^ = become re- 
dundant and only three energy-momentum conservation 
equations are required to be solved. 



account, fluid evolution depends only on seven compo- 
nents of the shear stress tensors. They are tt'^'^, tt"^^, 
^ry^ j^xx ^ ^yy ^ ^xy ^^^^ ^-q-q ^ Howcvcr, all the scvcn Com- 
ponents are not independent. Tracelessness, transver- 
sality to and the assumption of boost-invariance re- 
duces the independent components to three. Presently, 
we choose tt^^ tt^^ and tt^^ as the independent compo- 
nents. Relaxation equations for the independent compo- 
nents are given in the appendix [Cl (see Eas lC4IIC6ll . They 
are solved simultaneously with the three energy momen- 
tum conservation equations Egs lBlllBSl with inputs as 
discussed below. 




X (fm) X (fm) 



FIG. 1: (color online). Constant energy density contours in 
x-y plane at r=2.6 fm. The black lines are for ideal fluid 
{r]/s=0). The red and blue lines are for viscous fluid with 
ADS/CFT and perturbative estimate of viscosity, 77/5=0.08 
and 0.135. 

Heavy ion collisions are best described in [T,x,y^rj) 
coordinates, where r = ^/t"^ — z'^ is the longitudinal 
proper time and 77 = i In l^ti jg the space-time rapid- 
ity. r_L = {x, y) are the usual cartisan coordinate in 
the plane, transverse to the beam direction. Relevant 
equations concerning this coordinate transformations are 
given in the appendix 1X1 

Explicit equations for energy-momentum conservation 
in (r,x,y,77) coordinates are given in the appendix [Bl We 
note that unlike in ideal fluid, in viscous fluid dynam- 
ics, conservation equations (see Eas lBHlB3t contain addi- 
tional pressure gradients due to shear viscosity. Both T'^^ 
and T'^y components of energy-momentum tensor now 
evolve under additional pressure gradients. The right- 
most term of Eq lB3l also indicate that in viscous dynam- 
ics, longitudinal pressure is effectively reduced (note that 
the tt'''' component is negative). Since pressure can not 
be negative, shear viscosity is limited by the condition, 
p + r^T^m > 0. 

As evident from the Egs lBlllBSl in boost-invariant dis- 
sipative hydrodynamics, with shear viscosity taken into 



FIG. 2: (color online), same as Fig[l]but at time r=8.6 fm. 



IV. EQUATION OF STATE, VISCOSITY 
COEFFICIENT AND INITIAL CONDITIONS 

A. Equation of state 

One of the most important inputs of a hydrodynamic 
model is the equation of state. Through this input, the 
macroscopic hydrodynamic models make contact with 
the microscopic world. In the present demonstrative cal- 
culation we will show results for the QGP phase only. 
In the QGP phase, we use the simple equation of state, 
p — ^e, with energy density given as. 



£ = ^9,apT' (4.1) 

where gqgp = ggiuon + \gquark is the degeneracy factor for 
QGP. ggiuon — 2{helicity) x 8{color) is the degeneracy 
factor for gluons and gquark — 2{spin) x S{color) x 2{q + 
q) X Nf is the degeneracy factor for Nf flavored quarks. 
For Nf Ri 2.5, the degeneracy factor is gqgp = 42.25 
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B. Shear viscosity coefficient 

Shear viscosity coefficient (77) of QGP or sQGP is quite 
uncertain. In a strongly coupled QGP, shear viscosity 
can not be computed. Recently, using the ADS/CFT 
correspondence shear viscosity of a strongly coupled 
gauze theory, N=4 SUSY YM, has been evaluated, 77 = 
^N'^T^ and the entropy is given by s = ^N^T^. Thus 
in the strongly coupled field theory. 



contour plot of Vx at t:=8.6 fm 



(-) 



(4.2) 

ADS/CFT 47r 

Shear viscosity is quite uncertain in perturbative QCD 
also. At high temperature, shear viscosity, in leading log, 
can be written as [13, UHl , 

^3 

V = >^—. TT' (4-3) 

5 1115 

where g is the strong coupling constant. The leading 
log shear viscosity coefficient k depend on the number of 
fermion flavors {Nf). For example, for two flavored QGP, 

K = 86.47 and k = 106.7 for a three flavored QGP. With 

2 o 

entropy density of QGP, s = j^gggpT''. For two flavored 
QGP and ^0.5, the ratio of viscosity over the entropy, 
in the perturbative regime is estimated as, 



(-) 

Vs/j 



0.135, 



(4.4) 



' S / pert 

For lower as , perturbative estimation of rj/s could be 
even higher. 

Shear viscosity can also be expressed in terms of sound 
attenuation length, F^, defined as. 



sT 



(4.5) 



Ts is equivalent to mean free path and for a valid hy- 
drodynamic description F^/r << 1, i.e. mean free path 
is much less than the system size. Initial conditions of 
the fluid must be chosen carefully such that the validity 
condition Fs/r << 1 remains valid initially as well as at 
later time also. In the present work, we have treated vis- 
cosity as a parameter. To explore the effect of viscosity, 
we have used both the ADS/CFT estimate ?7/s=0.08 and 
perturbative estimate r;/s=0.135. We have also run the 
code with a higher value of viscosity 77/5=0. 2. 
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FIG. 3: (color online), contours of constant Vx in x-y plane 
at T—8.6 fm. The black lines are for ideal fluid (77/5=0). The 
red and blue lines are for viscous fluid with ADS/CFT and 
perturbative estimate of viscosity, 77/5=0. 08 and 0.135. 



contour plot of Vy at 1=8.6 fm 
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FIG. 4: (color online), contours of constant Vy in x-y plane 
at T=8.6 fm. The black lines are for ideal fluid (7;/s=0). The 
red and blue lines are for viscous fluid with ADS/CFT and 
perturbative estimate of viscosity, 7;/s=0.08 and 0.f35. 



C. Initial conditions 

Solution of Eqs lBHlB3l require initial conditions, the 
initial time t;, the transverse distribution of energy den- 
sity s{Xjjj) and the velocities Vx{x,y) and Vy(x,y). Fol- 
lowing initial transverse energy density is parame- 
terized geometrically. At an impact parameter &, trans- 
verse distribution of wounded nucleons Nwn{x, y, b) and 



of binary NN collisions Nbc{x, y, b) to are calculated in 
a Glauber model. A collision at impact parameter b is 
assumed to contain 25% hard scattering (proportional 
to number of binary collisions) and 75% soft scattering 
(proportional to number of wounded nucleons). Trans- 
verse energy density proflle at impact parameter b is then 
obtained as. 
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V. RESULTS 



£{x, y, h) = £0(0.75 X Nwn{x, y, b) + 0.25 x Nbc{x, y, b)) 

(4.6) 

with central energy density £o=30 GeVlfm"^. The equi- 
hbration time is chosen as Ti=0.6 fm [6]. The initial ve- 
locities Vx and Vy are assumed to be zero initially. 



VI. STABILITY OF NUMERICAL SOLUTIONS 



A. Evolution of the viscous QGP fluid 
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Temperature in x-x plane. y=0. 

ri/s=0 

ri/s=0.08 

ri/s=0 
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FIG. 5: (color online). Constant temperature contour in a; — 
T plane, for fixed y=0. The black, red and blue lines are 
for ideal, viscous fluid with r;/s=0.Q8 and viscous fluid with 
77/5=0.135. 



In dissipative hydrodynamics, one requires initial con- 
ditions for the viscous pressures also. Due to longitudinal 
boost invariance of the problem, we assume that viscous 
pressures have attained their boost-invariant values at 
the time of equilibration. Boost invariant values of the 
three independent shear stress-tensors can be easily ob- 
tained from Eqs [C4llC6l cr^^ ^ = 9 = ^ and a'^y 
( at the initial time r^, u^' = (1,0,0,0), Du^' = 0. The 
initial distribution of shear pressure tensors are then ob- 
tained as. 



7r-^(a;,y,6)= 2i^a-- = 2i^/n (4.7) 
T:yy (x, y, b) = 2'qayy = 2ri/n (4.8) 
7r="^(x,y,fe) = 2?7(t"^ = (4.9) 



As stated earlier, the viscous coefficient 77 is obtained 
using the relation, "q/s — const, const=0. 08, 0.135 or 0.2. 
For these values of shear viscosity, the validity condition 
Fg/r << 1 is satisfied initially. The validity condition is 
better satisfied at later time. 



The energy- momentum conservation equations IBlllBSl 
and the relaxation eg nations I C 2||C41 are solved simultane- 
ously using the code, AZHYDRO-KOLKATA, developed 
at the Cyclotron Centre, Kolkata. As mentioned earlier, 
we have solved the equations in the QGP phase only and 
did not consider any phase transition. In the following 
we will show the results for central Au+Au collisions (im- 
pact parameter 6 = fm). To understand the effect of 
shear viscosity, with the same initial conditions, we have 
solved the energy-momentum conservation equations for 
ideal fluid and viscous fluid. As mentioned earlier, we 
have considered two values of viscosity, the ADS/CFT 
motivated value 77/5=0. 08 and the perturbative estimate, 
77/5=0.135. 



Temperature in t-x plane. y=5 fm. 

Ti/s=0 

ri/s=0.08 

. 11/3=0.1 35 




FIG. 6: (color online), same as fig[5]but at y=5.0 

In Fig[Tl we have shown the contours of constant en- 
ergy density in x-y plane, after an evolution of 2.6 fm. 
The black lines are for ideal fluid evolution. The red and 
blue hnes are for viscous fluid with ADS/CFT (77/5=0.08) 
and perturbative (77/5=0.135) estimate of viscosity. Con- 
stant energy density contours, as depicted in FiglU indi- 
cate that with viscosity fluid cools slowly. Cooling gets 
slower as viscosity increases. Thus at any point in the x-y 
plane, energy density of viscous fluid is higher than that 
of an ideal fluid. At later time also, compared to an ideal 
fluid, viscous fluid evolve slowly. In FiglSj contours of 
constant energy density at time r=8.6 fm is shown. Here 
also we find than at any point energy density of viscous 
fluid is higher than its ideal counter part. The result is 
in accordance with our expectation. For dissipative fluid, 
equation of motion can be written as. 
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(6.1) 



Due to viscosity, evolution of energy density (or tem- 
perature) is slowed down. 

In Fig[3]and[4l we have shown the contour plot of the 
fluid velocity, and Vy, after evolution of 8.6 fm. As 
before the black lines are for the ideal fluid evolution. 
The red and blue lines are for viscous fluid with 7^/5=0.08 
and 0.135 respectively. Fluid velocities in viscous and 
ideal fluid differ very little. Even at late time, as shown 
in Figl3]and[3l we find that for 77/3=0. 08-0. 135, x and y 
component of the fluid velocity show marginal difference. 
However, there is an indication that in a viscous fluid, 
velocity grow faster than in ideal fluid. But as mentioned 
earlier, the difference is marginal. 
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FIG. 7: (color online). In the upper panel, temporal evolu- 
tion of the shear pressure tensor tt^^ at the fluid cell x=y=0 
is shown. In the lower panel, evolution of tt^^ at the fluid 
cell x=y=5 fm is shown. The black and red lines are for 
ADS/CFT motivated viscosity 77/5=0. 08 and perturbative es- 
timate 77/5=0.135 respectively. 

As seen in FiglTJH in viscous dynamics, QGP fluid 
evolves slowly. Thus life-time of the QGP phase is en- 
hanced in viscous dynamics. To obtain an idea about the 
enhanced life-time, in Fig[51 we have shown the constant 
temperature contours m t — x plane , at a fixed value 
of y=0 fm. As seen in Fig|5l temperature evolves slowly 
in a viscous fluid and life-time of the QGP phase is ex- 
tended. For small viscosity 77/5=0. 08-0. 135, the increase 
is not large. At the center of the fluid, for 77/5=0. 135, 
QGP life-time is increased approximately by 5% only. It 
is even less for the ADS/CFT estimate of viscosity. How- 
ever, enhancement of QGP life-time depends on the fluid 



cell position. It could be more. In Figl6l constant tem- 
perature contours at y=5 fm is shown. For 77/5=0.135, at 
x=0,y=5 fm, the QGP life-time is enhanced by ~ 10%. 
We conclude that in a viscous dynamics, with moder- 
ate viscosity 77/5=0.08-0.135, QGP life-time could be en- 
hanced by 5-10%. Enhanced lifetime of QGP in a viscous 
fluid can have significant effect on observables produced 
early in the collisions e.g. direct photon production or in 
J/t/; suppression. 




X (fm) 



FIG. 8: (color online). In panel (a) and (b), contours of 
constant pressure tensor tt^^ at initial time Ti— 0.6 fm and at 
time r=2.6 fm is shown. In panel (c) and (d) same results for 
shear pressure tensor tt^" is shown. 



B. Evolution of shear pressure tensors 

We have assumed that initially the shear pressure ten- 
sors TT^^, w'^y and ir^y attained their longitudinal boost- 
invariant values. As the fluid evolve, pressure tensors also 
evolve. Here we investigate the evolution of shear pres- 
sure tensors with time. In the top panel of Fig[7] evolu- 
tion of shear pressure tensor tt^^ at the fluid cell position 
x=y=0 is shown. The black hne is for the ADS/CFT 
motivated viscosity, 77/5=0.08 and the red line is for the 
perturbative estimate of viscosity 7//5=0.135. Just after 
the start of the evolution the shear pressure tensor tt^^ 
increases, but for a short duration and then steadily de- 
creases with time. By 4 fm of evolution, tt^^ at the center 
of the fluid reduces to negligibly small values. Identical 
behavior is seen for the shear pressure tensor tt^^. In 
the bottom panel of Fig[7] we have shown the evolution 
of the third independent shear pressure tensor tt^^. Ini- 
tially TT^y is zero. As the fluid evolve, it grow in the 
negative direction. We find that at the centre of the fluid 
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(x=y=0), it never grows. In Fig[7l temporal evolution 
of TT^*' at the fluid cell position x ~ y — 5fm is shown. 
From the initial zero value, tt^^ rapidly increases in the 
negative direction. It reaches its maximum around t «1 
fm and then decreases again. We also note that tt^** never 
grows to large values. Compared to tt^^ or tt^^ stress ten- 
sor TT^y is negligible. The results indicate that in a QGP 
fluid, viscous effect persist for a short duration (3-4 fm) 
only. At late time the fluid evolve essentially as an ideal 
fluid. The result is understandable. Shear viscosity de- 
pend strongly on temperature (77 oc T"^). As the fluid 
cools, effect of viscosity decreases rapidly. 



C. Entropy generation 

In a viscous fluid dynamics, entropy is generated. We 
can easily calculate the entropy generated during the evo- 
lution. 



2r]T 



(6.2) 
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FIG. 9: Evolution of average entropy with time, for two val- 
ues of viscosity, the ADS/CFT motivated viscosity 77/5=0. 08 
and perturbative estimate 77/3=0. 135 are shown. 



To show the spatial distribution of the stress tensors, 
in FiglH TT^^ and tt^^ at initial time Ti=0.6 fm and after 
an evolution of r = 2.6 fm are shown. As shown earlier, 
TT^^ and also tt^^ rapidly decreases with time. By 2 fm 
of evolution they are reduced by approximately by a fac- 
tor 6. It is also interesting to note that the initial x-y 
symmetric distribution of tt^^ and tt^^ quickly evolves to 
asymmetric distribution. With time tt^^ evolves faster in 
the x-direction than in y-direction. Similarly, tt^^ evolve 
faster in the y-direction than in the x-direction. For cen- 
tral collisions the asymmetric evolution of tt^^ and tt^^ 
counter balance each other. As shown in Fig[T] and [H 
the contour plots of energy density do not show any in- 
dication of asymmetry even at late time. However, the 
asymmetric pressure tensors can have important effects 
on elliptic flow of observables produced early in the col- 
lisions, say in elliptic flow of direct photons. 



Evolution of spatially averaged entropy is shown in 
FigEl for the two values of viscosity coefflcients 77/5=0.08 
and 0.135. As expected, entropy generation is more if vis- 
cosity is more. For both the values of viscosity, we find 
that entropy generation saturates after « 3 fm of evolu- 
tion. It is expected also. As shown previously, viscous 
fluxes reduces to very small values after r=3 fm. Natu- 
rally, entropy generation is negligible thereafter. 



Temperature in x-x plane, y= 5m. 

ideal fluid 

viscous fluid (1st order) 

viscous fluid (2nd order) 



E 
> 




X (fm) 



FIG. 10: (color online) constant temperature contours in 
X — T plane at y=5 fm. The black lines are for ideal fluid. 
The red and blue lines are for viscous fluid in 1st order and 
in 2nd order theory respectively. 77/3=0. 135. 



D. 1st order theory vs. 2nd order theory 

As mentioned earlier, 1st order theory of dissipative 
hydrodynamics is acausal, signal can travel faster than 
light. This is corrected in 2nd order theory, but we 
have to pay the price, relaxation equations for dissipa- 
tive fluxes are required to be solved. It is interesting to 
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FIG. 11: Evolution of average entropy production in a 1st 
order (solid line) and 2nd order (dashed line) theory. 2nd 
order theory generate more entropy. 



compare the difference we can expect in a first order the- 
ory and in a 2nd order theory of dissipation. In Fig. 1101 
we have shown the contours of constant temperature in 
a;— r, for a fixed y = 5fm. The black hues are for an ideal 
fluid. The red hues are for a viscous fluid treated in the 
1st order theory. The blue lines are for viscous fluid in 
2nd order theory. In 2nd order theory fluid evolve more 
slowly than in a flrst order theory. Entropy generation is 
also more in a 2nd order theory. In FigfTTl average en- 
tropy evolution with proper time is shown, both for the 
1st order theory (the solid line) and the 2nd order theory. 
In 2nd order theory, approximately 80% more entropy is 
generated. 



VII. TRANSVERSE MOMENTUM AND 
ELLIPTIC FLOW OF QUARKS 

Presently we can not compare predictions from viscous 
hydrodynamics with experimental data. Hadrons are not 
included in the model. The initial QGP fluid evolve and 
cools but remain in the QGP phase, it did not undergo 
a phase transition to hadronic gas. However, from the 
momentum distribution of quarks we can get some idea 
about the viscous effect on particle production. Viscos- 
ity generates entropy, which will be reflected in enhanced 
multiplicity. We use the standard Cooper-Frey prescrip- 
tion to obtained the transverse momentum distribution 
of quarks. In Cooper-Frey prescription, particle distribu- 
tion is obtained by convoluting the one body distribution 
function over the freeze-out surface, 



(7.1) 



where dE^ is the freeze-out hyper-surface and f{x,p) is 
the one-body distribution function. Now in a viscous dy- 
namics, the fluid is not in equilibrium and f{x,p) can not 
be approximated by the equilibrium distribution func- 
tion, 



1 



expl 



(7.2) 



with inverse temperature (3 = 1/T and chemical poten- 
tial II. In a highly non-equilibrium system, distribution 
function f{x,p) is unknown. If the system is slightly 
off-equilibrium, then it is possible to calculate correction 
to equilibrium distribution function due to (small) non- 
equilibrium effects. Slightly off-equilibrium distribution 
function can be approximated as. 



fix,p)^ f^°\x,p)[l + ^ix,p)l 



(7.3) 



(j){x,p) is the deviation from equilibrium distribution 
function /f"^. With shear viscosity as the only dissi- 
pative forces, cf>(x,p) can be locally approximated by a 
quadratic function of 4-momentum, 



Without any loss of generality 



2{e + p)T^ 



(7.4) 

can be written as 
(7.5) 



completely specifying the non-equilibrium distribution 
function. As expected, correction factor increases with 
increasing viscosity. We also note that non-equilibrium 
correction is more on large momentum particles. The 
effect of viscosity is more on large momentum parti- 
cles. The correction factor reduces if freeze-out occurs 
at higher temperature. 

With the corrected distribution function, we can cal- 
culate the quark momentum spectra at freeze-out sur- 
face E^. In appendix iDl relevant equations are given. 
The quark momentum distribution has two parts, (i) 
dyd^PT ' obtained by convoluting the equilibrium distribu- 
tion function over the freeze-out surface and (ii) ^^'ij^ , 
obtained by convoluting the correction to the equilibrium 
distribution function over the freeze-out surface. Since 
the correction factor is obtained under the assumption 
that non-equilibrium effects are small, <j){x,p) « 1, it 
necessarily imply that, « The ratio. 



dN<"i 



dyd'-'pT 

dN'"i ' 
dyd'-^PT 



(7.6) 



could at best be unity or less. If the ratio exceeds 
unity, it will imply that non-equilibrium effects are large 
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FIG. 12: Ratio of quark spectra with non-equilibrium distri- 
bution function to that with equilibrium distribution function. 



and the distribution function f{x,p) can not be approx- 
imated as in EqO Using AZHYDRO-KOLKATA, we 

and 



have simulated a b=6.5 fm Au-|-Au collision. 



dN" 



dN"" 
dyd^PT 



-3--J2 — at freeze-out temperature Tf=160 MeV are cal- 

dyd'^PT 

culated. The ratio for r//s=0.08,0.135 and 0.2, are 

shown in FigdU With ADS/CFT estimate of viscosity, 
77/5=0.08, non-equilibrium correction to particle produc- 
tion become comparable to equilibrium contribution only 
beyond pt—5 GeV. However, with perturbative estimate, 
?7/s=0.135, non-equilibrium correction become compara- 
ble to or exceeds the equilibrium contribution at px ~ 
4 GeV. Pt range is further reduced for higher viscosity 
r]/s=0.2. Thus with perturbative estimate of viscosity 
(77/s = 0.135 — 0.2), hydrodynamic description remain 
valid upto transverse momentum pT ~ 3.5-4 GeV. 

In Fig |131 we have compared the transverse momen- 
tum spectra of quarks in ideal hydrodynamics with that 
in a viscous dynamics. In Fig[T31 the dotted line is the 
spectra obtained in ideal dynamics {rj/s — 0). The pT 
spectra in viscous dynamics are shown by black lines. 
We have shown the spectra for three values of viscos- 
ity r//s=0. 08, 0.135 and 0.2. Compared to ideal dynam- 
ics, quarks yield in viscous dynamics increases. The in- 
crease is more at large pt- For low values of viscosity 
the increase is modest, a factor of 2 at pt—3 GeV. But 
yield increase by a factor or 4(10) if viscosity increases 
to 77/5=0.135 (2). Please note that even though we have 
shown Pt spectra upto 5 GeV, for 77/5=0.2 and 0.135, 
hydrodynamic description fails beyond pT ^ 3.5 and 4 
GeV. 

We have also studied the effect of viscosity on quark 
elliptic flow. Effect of viscosity is very prominent on 
elliptic flow. In Fig[T31 pT dependence of elliptic flow 
of quarks, in a b=6.5 fm collision is shown. The black 
line is V2 in ideal dynamics. In ideal dynamics, ellip- 
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FIG. 13: Qurak transverse momentum spectra at freeze-out 
temperature of 160 MeV. The dotted line is the quarks spectra 
in ideal hydrodynamics. The solid lines (top to bottom) are 
in viscous dynamics with 77/3=0. 08, 0.135 and 0.2. 



tic flow continually increases with pT- It well known, in 
contrast to experiments, where elliptic flow saturates at 
large pT, in ideal hydrodynamics, elliptic flow continue 
to increase with pT- Indeed, this is a major problem in 
ideal hydrodynamics. The renewed the interest in dis- 
sipative hydrodynamics is partly due to the inability of 
ideal hydrodynamics to predict the trend of elliptic flow 
in Au-|-Au collisions. In FigfT4l the blue lines are V2 
in viscous dynamics with 77/5=0. 08, 0.135 and 0.2 respec- 
tively. In a viscous dynamics, pt dependence of V2 is 
drastically changed. In contrast to ideal dynamics where 
V2 continue to increase with pT, in viscous dynamics, V2 
continue to increase only upto pt ~ 1.5 — 2GeV. There- 
after V2 decreases. For perturbative estimate of viscos- 
ity 77/5=0.135 and beyond, V2 even become negative at 
large pt ■ Veering about of V2 after pT ~ 1.5-2 GeV is 
due to viscous effect only or more explicitly due to the 
non-equilibrium correction to the equilibrium distribu- 
tion function. This is clearly manifested from the red 
lines in FigHD The red lines are calculated ignoring the 
non-equilibrium corrections to the equilibrium distribu- 
tion function. If non-equilibrium correction is ignored, in 
viscous dynamics also, V2 continue to increase with pT, 
albeit its magnitude is reduced compared to ideal dy- 
namics. The result is very important. It imply that the 
experimental trend of elliptic flow (saturation at large 
Pt) could only be explained if the QGP fluid is viscous. 
An ideal QGP, will not be able to explain the saturation 
trend of the experimental data. 

As stated earlier, non-equilibrium correction to the 
equilibrium distribution function depends on the freeze- 
out condition. To show the effect of freeze-out condition, 
on V2, in Fig |15l we have shown V2 for a values of freeze- 
out temperature Tp =160,150,140,130 and 120 MeV. As 
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FIG. 14: (color online) Elliptic flow as a function of transverse 
momentum. The black line is V2 in ideal hydrodynamics. The 
blue lines are 112 in viscous dynamics with viscosity to entropy 
ratio 77/5=0.08,0.135 and 0.2 (top to bottom) respectively, 
including the correction to equilibrium distribution function. 
The red lines are same as the blues lines but ignoring the 
non-equilibrium correction to the distribution function. 

freeze-out occur at higher and higher temperature, the 
veering of V2 takes place at larger and larger px and for 
Tf=120 MeV, the eUiptic flow saturates. The result is 
understood easily. With decreasing freeze-out tempera- 
ture, the fluid evolves for longer time, the shear stress- 
tensor's at the freeze-out surface is reduced and the non- 
equilibrium correction, proportional to shear stress ten- 
sors, decreases. 



VIII. STABILITY OF NUMERICAL 
SOLUTIONS IN AZHYDRO-KOLKATA 

Before we summarise our results, we would like to 
comment on the stability of numerical solutions in 
AZHYDRO-KOLKATA. As indicated above, with shear 
viscosity as the only dissipative force, boost-invariant 
causal hydrodynamics require simultaneous solution of 
six partial differential equations. Numerical solution of 
six partial differential equations is non-trivial and it is 
important to check for the numerical stability of the 
solutions. Analytical solutions of viscous hydrodynam- 
ics, even in restrictive conditions are not available, and 
we can not check the solutions against analytical re- 
sults. However, we can check for the stability of the 
numerical solutions. The standard procedure of check- 
ing the numerical stability is to change the integration 
step lengths and look for the difference in the solution. 
In Fig |161 for viscosity 77/5=0. 135, we have shown the 
constant temperature contours in a; — t plane at a fixed 
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FIG. 15: Dependence of elliptic flow on the freeze-out tem- 
perature. The solid lines (from bottom to top) are elliptic 
flow («2) in viscous dynamics with rF = 160,150,140,130 and 
120 MeV respectively. 



2 



Au+Au@b 


=0 


dx=dy=0.2,dT=0.02 
HY=H\/=n 1 Ha-=n ni 




"Imp 
















16 






















■0t3^ 






I 





I ' I" — ' T — — r ' — I ' 1 

2 4 6 8 10 

X (fm) 



FIG. 16: (color online) constant temperature contours in 
X — T plane at a fixed y=0 fm. The black lines are obtained 
with integration step lengths dx=dy=0.2 fm and dr=0.02 fm. 
The blue lines are obtained with integration step lengths, 
dx=dy=0.1 fm and dr=0.01 fm. Halving the step lengths do 
not change the evolution. The numerical solutions are stable. 



y=0 fm. The black and blue lines are obtained when 
integration step lengths are dx=dy=0.2fm,c?T=0.02 fm, 
and dx=dy=0.1 fm,c?T=0.01fm respectively. Evolution 
of QGP fluid donot alter by changing the step lengths, 
the solutions are stable against mesh size. 
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IX. SUMMARY AND CONCLUSIONS 

In Israel-Stewart's 2nd order theory of dissipative rel- 
ativistic hydrodynamics, wc have studied evolution QGP 
fluid. In 2nd order theory, in addition to usual thermo- 
dynamic quantities e.g. energy density, pressure, hydro- 
dynamic velocities, dissipative flows arc treated as ex- 
tended thermodynamic variables. Relaxation equations 
for dissipative flows are solved, simultaneously with the 
energy-momentum conservation equations. This greatly 
enhances the complexity of the problem. Altogether 14 
partial differential equations arc required to be solved. 
We simplify the problem to solution of six partial differ- 
ential equations by considering the evolution of baryon 
free QGP fluid with longitudinal boost-invariance. We 
also consider dissipation due to shear viscosity only, dis- 
regarding the bulk viscosity and the heat conduction 
(for a baryon free QGP fluid they do not contribute). 
The six partial differential equations are solved using 
the code AZHYDRO-KOLKATA, developed at the Cy- 
clotron Centre, Kolkata. 

To bring out the effect of viscosity, we have considered 
the evolution of ideal as well as viscous QGP fluid. Both 
ideal and viscous fluid are initialized similarly, at initial 
time Ti=0.6 fm, the central entropy density is 110 fm~^. 
Viscous dynamics require initial conditions for the shear- 
stress tensor components. It is assumed that at the equi- 
libration time, the shear stress tensors components have 
attained their boost-invariant values. 

Explicit simulation of ideal and viscous fluids confirms 
that energy density of a viscous fluid, evolve slowly than 
its ideal counterpart. Thus in a viscous fluid, lifetime of 
the QGP phase will be enhanced. Transverse expansion 
is also more in viscous dynamics. For a similar freeze-out 
condition freeze-out surface is extended in viscous fluid. 

As the fluid evolve, shear pressure tensors also evolve. 
Explicit simulations indicate that shear pressure tensors 
TT^^ and TT^^ which are initially non zero, rapidly de- 
creases as the fluid evolve. By 3-4 fm of evolution they 
reduced to very small values. The other independent 
shear tensor tt^^ is zero initially. At later time it grow 
in the negative direction but never grow to large value 
and is always order of magnitude smaller than the stress 
tensors (tt"^^ and tt^^). Spatial distribution of shear pres- 
sure tensors tt^^ and tt^^ reveal an interesting feature 
of viscous dynamics. Initially tt^^ and tt^^ have sym- 
metric distribution. As the fluid evolve, pressure tensors 
quickly become asymmetric, e.g. tt^^ evolve faster in the 
x-direction than in the y-direction, tt^^ evolve faster in 
y direction than in x-direction. However, in a central 
collision, wc; did not see any effect of asymmetry in the 
energy density distribution. In a central b=0 collision, 
the two opposite asymmetry cancels each other. 

We could not study effect of shear viscosity on par- 
ticle production. However, we have explored the effect 
of viscosity on parton momentum distribution and ellip- 
tic flow. We have simulated b=6.5 fm Au-|-Au collision. 
Using the Cooper-Frey prescription, transverse momen- 



tum spectra as well as elliptic flow of quarks at freeze-out 
temperature of Tp^WO MeV are obtained. Viscous dy- 
namics flattens the quark yield at large px- At pt=3 
GeV, even a small viscosity, t]/s=0.8, increase the yield 
by a factor of 2. The increase is even more if viscosity is 
large. Viscous effect is most prominent on elliptic flow. 
In ideal hydrodynamics, elliptic flow continue to increase 
with pt- But in viscous dynamics V2 veer about around 
Pt=1-5-2 and even become negative at large pT- With 
appropriate choice of viscosity, freeze-out condition, el- 
liptic flow show saturation. The saturation effect is es- 
sentially due to non-equilibrium correction to the equi- 
librium distribution function and can not be mimicked 
in an ideal hydrodynamics. Only in viscous dynamics, 
saturation of elliptic flow can be explained. 



APPENDIX A: COORDINATE 
TRANSFORMATIONS 

Instead of Cartesian coordinates = {t, x, y, z) we 
use curvilinear coordinates in longitudinal proper time 
and rapidity, = (r, x, y, t]): 



t = Tcoshr/; 
z = Tsinhr/; 

The differentials 



1, t^z 

- In 

2 t-z 



dt = dr cosh i] + d-qr sinh r], 
dz = dr sinh rj + drjT cosh rj, 

and the metric tensor is easily read off from 

ds'^ = g^^dxf^dx" = dt^ - dx"^ - dy"^ - dz"^ 



(Al) 
(A2) 



(A3) 
(A4) 



Qmndx'^dx'^ = dr^ - dx^ - - r2d??^(A5) 



namely 



9mn 



10 
0-10 
0-10 
-r^ 



/ 1 
0-10 
0-1 
VO -1/t2 
(A6) 

In curvilinear coordinates we must replace the partial 
derivatives with respect to a;'' by covariant derivatives 
(denoted by a semicolon) with respect to x^: 

- fin^ik _ _ 

rpik _ |_ pi rpmk i rpim-pk 

dxP ^ 

The only non- vanishing Christoffel symbols are 

r>^ = r;j, - i/r. (A?) 



The hydrodynamic 4-velocity = j{l,Vx,Vy,Vz) is 
transformed to u™ = 7(1, Uj,, %, 0), with 7j_ = 
l/-y/l— From here on, we drop the bars over tensor 
components in x-coordinates for simplicity. 
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The projector can be easily calculated, 



(\ 



V 



-7i 

-l\Vy 
















^1 



(A8) 



In {T,x,y,ri) coordinate system, the convective time 
derivative can be obtained as, 



D = u - d = 7(9r + Vxdx + Vydy). 



(A9) 



For future reference,we also write down the the scalar 
expansion rate 



d-u = dru^ + + dyu^ + — (AlO) 



APPENDIX B: ENERGY-MOMENTUM 
CONSERVATION 

With longitudinal boost-invariance the energy- 
momentum conservation equations T"^"'-n = yield 



drT^^ + d,{T^^v,) + dyiT^^Vy) = -{p + t^tt"") 

drf-y + d.,{f^yv,) + dyif-yvy) = ~d,{n-y - ^yy,) - dy{p + %yy - ^^yvy) 



(Bl) 
(B2) 
(B3) 



where i"" = rA"", p = rp, and = T^^/T^^, Vy = 

The components of the energy momentum tensors, in- 
cluding the shear pressure tensor are, 

r""= (e+j5)7i -P + tt"" (B4) 
T-^'- (e+p)7i^'. + ^"" (B5) 
ry ^ {^e + p)^lvy + (B6) 

In causal dissipative hydrodynamics, energy momen- 
tum conservation equations are solved simultaneously 
with the relaxation equations. Given an equation of 
state, if energy density (e) and fluid velocity [vx and 
Vy) distributions, at any time r.^ are known, Eqs lBllB2l 
and IB3I can be integrated to obtain e, Vx and Vy at the 
next time step Ti+i. While for ideal hydrodynamics, this 
procedure works perfectly, viscous hydrodynamics poses 
a problem that shear stress-tensor components contains 
time derivatives, 9r7±, drU^, drU^ etc. Thus at time step 
Ti one needs the still unknown time derivatives. Numer- 
ically, time derivatives at step Ti could be obtained if ve- 
locities at time step and r^+i are known. One possible 
way to circumvent the problem, is to use time derivatives 
of the previous step, i.e. use velocities at time step Ti-i 
and Ti to calculate the derivatives at time step [l8l|. 
The underlying assumption that fluid velocity changes 
slowly with time. In 1st order theories, this problem is 
circumvented by calculating the time derivatives from the 
ideal equation of motion , 



e+p 

De = -(e+p)V^M^. (B8) 



With the help of these two equations all the time 
derivatives can be expressed entirely in terms of spa- 
tial gradients [l^, . 1st order theories are restricted 
to contain terms at most linear in dissipative quantities. 
Neglect of viscous terms can contribute only in 2nd or- 
der corrections, which are neglected in 1st order theories. 
While the procedure is not correct in 2nd order theory, 
we still use it in the present calculations. The alternative 
procedure of using the derivative of earlier time step is 
not correct either. 



APPENDIX C: RELAXATION EQUATIONS FOR 
THE VISCOUS PRESSURE TENSOR 

Being symmetric and traceless, the viscous pressure 
tensor tt'"' has 9 independent components. The as- 
sumption of boost invariance reduces this number by 
3 (V^™m''^=0, m 7^ 7]). The transversality condition 
UtoTt'"" = eliminates another three components ( 
vanish and thus yield no constraint). Thus, with boost- 
invariance the viscous pressure tensor has only three in- 
dependent components. As seen in Eqs lBllB2l and [B3l in 
a boost-invariant evolution only seven pressure tensors 

^TT^ j^xx ^ ^yy ^ ^rjr) ^ j^rx ^ ^ry ^xy ^^^.^ importance. 

Only three of these seven are independent. In an ear- 
lier publication [T^, we have debated about the choice 
of the independent components and suggested use of ei- 
ther (tt^^, tt'''', a = vr^^ - ttW) or (7r^^,7rW,7r^^, tt^w) 
(which will require solution of an additional relaxation 
equation) as choice of independent components. How- 
ever, while computing we find that the three pressure 
tensors tt^^ and tt^^ and tt^^ as independent components 
are computationally more convenient. The choice has 
the advantage that the dependent shear stress tensors 
can be obtained from the 3 independent stress tensors by 
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multiplying them by fluid velocity, Vx and Vy (see Eqs. 
ICTHClOp . In any other choice of independent components 
(e.g. 7r'^'^,7r'''',A = tt^^ — tt^^), the evaluation of depen- 
dent stress tensors requires division by fluid velocities. 
Since initially, fluid velocities are assumed to be zero and 
they grow slowly, these choices will involve division by 
very small numbers. Unless proper care is not taken, di- 



vision by small numbers can lead to unrealistically large 
values for the dependent stress tensors and ruin the com- 
putation. 

The relaxation equations for the independent shear 
stress tensors tt^^, tt^^ and tt^**, in (r,x,y,77) co-ordinate 
can be written as, 



-— (7r^^-2r7CT='^) 
-—(n-'y-2na-'y) 



(CI) 
(C2) 
(C3) 



r 



where is the relaxation time, T,r = 277/32 (see 
Ea l2.2ip . In ultra-relativistic limit, for a Boltzman gas, 
P2 can be evaluated, /32 ~ ^ where p is the pressure 
[T3| . In the present paper, we use this limit to obtain the 
relaxation time r^. 

The viscous pressure tensor relaxes on a time scale 
to 277 times the shear tensor a'^'^ = V'^^u'^\ The xx, 
yy and xy components of the shear tensor cr^" can be 
written as 



rVV 



-dyuy - u^Duy 



/svvo 



(C4) 
(C5) 



A'^ye (C6) 



The dependent shear stress tensors can easily be ob- 
tained from the independent ones as. 



(C7) 
(C8) 
(C9) 

^l)nyy 



TT = 


XX 

VxTT 


+ VyTT^'y 


TT •' — 


vx^^y 


H- vyTTyy 


^'^'^ 
TT — 


vItt-^ 


+ vl^yy H 


T TT ' ' = 


-(1- 





T-xy 

XX -{l-vl 

+ 2VxVyTT^y 



(CIO) 



The expressions for the convective time derivative D 
and expansion scalar 9 — dii, in (T,x,y,77) are given in 
Eqs. Egiand lUOl 



APPENDIX D: PARTICLE SPECTRA 

With the non-equilibrium distribution function thus 
specified, it can be used to calculate the particle spectra 



from the freeze-out surface. In the standard Cooper- Frye 
prescription, particle distribution is obtained as. 



E 



dN dN 



d^p dyd^pT Jt. 



d^^p^f{x,p) (Dl) 



In (T,x,y,ris) coordinate, the freeze-out surface is pa- 
rameterised as, 



E'" = (t/(x, y) cosh 77s, a;,y,T/(a;,y)sinhr/s), (D2) 
and the normal vector on the hyper surface is, 



dTf dTf 

dT.^, = (cosh f]s , - - — , - - — , - sinh rjs )Tfdxdydris 
dxf dyf 

At the fluid position (T,x,y,ris) the particle 4- 
momenta are parameterised as, 



pf^ = [mTCOsh{ris — Y) , p'^ , py , mTsinh{ris — Y)) (D4) 
The volume element p^dH^^ become, 

p'^dT.^ — (rriTCOshir] — Y) ~ pT-'^TTf)Tfdxdydr] (D5) 

Equilibrium distribution function involve the term 
which can be evaluated as. 



P^Ufj, j(rnTCosh(ri — Y) ~ vt-Pt ~ fJ-fj) 



(D6) 



T T 

The non-equilibrium distribution function require the 
sum pf^p^TTf^iy, 

Pfj.Pu'n''^'^ — aiCOsh?{r] — Y) + a2Cosh{rj — Y)+a^ (D7) 
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with 



grating over spatial rapidity one obtains, 



a2 = -2mTiPx'^^'' +PyTT'''") 



(D8) 
(D9) 



Inserting all the relevant formulas in Eq IDll and inte- 



dN 



dyd^PT dyd^pT dyd^px 



with, 



(DU) 



dN""! 
dyd^pj 



dxdyTf[mTKi{n(3) — pT'^TTfKQ{nf3)] 



-pT.VTT/{yif2M) + aaifiM) + (y + a3)Ko{nP)}] 



r 



(D12) 



(D13) 



where Kq, Ki, K2 and if 3 are the modified Bessel func- 
tions. 

We will also show results for elliptic flow V2- It is de- 
fined as. 



V2 



f27r dN 
J dyd'^pT 



cos{2(f))d(f) 



/c 



27r dN 



(D14) 



JO dyd^pT ^ 

Expanding to the 1st order, elliptic flow as a function 



of transverse momentum can be obtained as, 



d^N" 



V2{Pt) = V2^{pt) 1 - 



Pt dpTd4> 



J d(j)cos{2(j)) 



PTdpTd4> 



I' 



Pt dpT dq 



Pt dpT d(j) 

(D15) 

where v'^'^ is the elliptic flow calculated with the equi- 
librium distribution f^"^ . 
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